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The Magnetospheric Multiscale (MMS) mission consists of four identically in- 
strumented, spin-stabilized observatories elliptically orbiting the Earth in a tetra- 
hedron formation. A requirement for the operational success of the mission is 
the ability for the on-board systems to deliver precise maneuver adjustments. A 
six degree-of-freedom (6-DOF), closed-loop control system was developed that 
tracks a time-varying, inertial velocity-target with less than 1% error down to a five 
millimeter-per-second lower-threshold (3a). This level of performance is achieved 
in-part through integrated and dynamically-compensated accelerometer feedback 
with micro-gravity resolution. System performance is bounded through an exten- 
sive Monte Carlo simulation campaign that exercises the multi-body dynamics and 
non-linear sensitivities, and supported by some initial flight-results. 


INTRODUCTION 

The Magnetospheric Multiscale (MMS) mission, launched on March 13, 2015, is the fourth mis- 
sion of NASA’s Solar Terrestrial Probe program. The MMS mission consists of four identically 
instrumented observatories that function as a constellation to provide the first definitive study of 
magnetic reconnection in space. 

Since it is frequently desirable to isolate electric and magnetic field sensors from stray effects 
caused by the spacecraft’s core -body, the suite of instruments on MMS includes six radial and 
two axial instrument-booms with deployed lengths ranging from 5-60 meters (see Figure 1). The 
observatory is spin-stabilized about its positive z-axis with a nominal rate slightly above 3 rev/min 
(RPM). The spin is also used to maintain tension in the four radial wire -booms. 

Each observatory’s Attitude Control System (ACS) consists of digital sun sensors, star cameras, 
accelerometers, and mono-propellant hydrazine thrusters — responsible for orbital adjustments, at- 
titude control, and spin adjustments. The sections that follow describe performance requirements, 
the hardware and algorithms used for 6-DOF estimation, and then similarly for 6-DOF control. The 
paper concludes with maneuver performance based on both simulated and on-orbit telemetry. 

Performance Requirements 

Once in science mission orbits, the four 0.12-km diameter observatories plan to form a tetrahe- 
dron with as little as 4-km of separation between spacecraft. The need to maintain a sufficiently 
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Figure 1. MMS Observatory (Stowed / Fully-Deployed) 


accurate spatial and temporal formation must be balanced against the logistical constraints of ex- 
ecuting overly-frequent maneuvers on a fleet of spacecraft. These two considerations make for an 
extremely challenging orbit design problem. 1,2 


With a stated operational goal of maneuvering the fleet no 
more often than once every two weeks (on average), MMS 
flight dynamics specialists derived the formation maintenance 
maneuver execution requirements shown to the right. 


ATTITUDE, RATE, AND ACCELERATION 
Star Sensor 

The //A SC Star Sensor System 
(STS), provided by the Technical 
University of Denmark, consists 
of internally redundant electronics 
housed within a single enclosure 
that interfaces with four Charge- 
Coupled Device camera head units. 

The STS provides time-stamped at- 
titude quaternion data packets at a 
4 Hz telemetry rate. It has a 3er full 
performance transverse and bore- 
sight axis accuracy of 60 arcsec and 
200 arcsec, respectively. The STS 
has a spin rate capability of up to 
4 RPM. Figure 2 shows a 752x580 
pixel image taken April 27, 2015 by 
Observatory 3 (MMS-3) using STS camera head unit B. 
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Figure 2. DTU /jASC Image (OBS 3, CHU B) 


Multiplicative Extended Kalman Filter (MEKF) 

Attitude and angular rate of the spacecraft with respect to the ECI J2000 inertially fixed reference 
frame are determined on-board by processing the four quaternion solutions produced by each star 
sensor camera-head unit (CHU) at 4 Hz. The CHU measurements are transformed into the com- 
mon spacecraft body-fixed frame and then combined in a Multiplicative Extended Kalman Filter 
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(MEKF). Due to limitations in the processing-power of the flight computer (Motorola RH-CF5208 
Coldfire), neither the sun sensor nor the acceleration measurements arc included in the filter’s com- 
putations. Also, due to an early programmatic design-trade, MMS does not include gyroscopes 
in its complement of sensor. As a result, the MEKF does not use the standard gyro-replacement 
model and instead follows the method of Thienel, Mark Icy and Harman.’ A brief overview of the 
MMS on-hoard algorithm follows. The presentation assumes familiarity with the application of 
Kalman filters for spacecraft attitude estimation. As a result, only the elements particular to our 
system model will be discussed. A more complete explanation of the MEKF mechanics (covari- 
ance/gains/update equations, etc.) can be found in the references. 4 ' 5 

States The MMS filter uses a seven-component vector consisting of the body’s attitude quater- 
nion q and the angular rate uo for its full-state X(q, uo), and a reduced six-component based on the 
Gibbs vector a attitude representation for its error-state x(a, Suo). The relationship between the 
error-quaternion Sq (defined as the difference in attitude between the true attitude q at the time of 
measurement and the propagated estimate q _ ) and the Gibbs vector is, to lst-order, 

a' 

2 ( 1 ) 

1 

where the operator ® denotes the quaternion product. 5 

Measurements for Update The star sensor measurement equation is simply the measured quater- 
nion of the individual camera-heads transformed to the base-body frame, 

y = h(q, t) + noise = q chu (2) 


Sq = q<g) (q ) 1 « 


Since the MEKF error states are “reset” after an update and completely absorbed in the full-state 
estimate (i.e. a + = 0), the measurement residuals p is equal to the quaternion product between the 
a priori attitude quaternion estimate q and the current measurement. 


P = y - h (q) - H a d 4 (3) 

= Q ch u = 2 • dq (1:3) • sign (Sq 4 ) (4) 

where q 4 is the quaternion scalar, and H a is a component of the linearized measurement sensitivity 
matrix, 
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and I is a 3 x 3 identity matrix. 

Model The full non-linear state dynamics, 
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arc the quaternion kinematics and Euler’s rotational equation , driven by the control-input torque 
r(u) and white noise w. The superscript expression [] T is the matrix transpose, and [] x denotes a 
skew-symmetric 3x3 matrix formed from the base vector. For example, 
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Also, note that the presence of the body’s second mass-moment of inertia I in Eq. (7) requires 
accurate knowledge of the inertia mass properties for rate estimation — something that is typically 
avoided via direct substitution of angular rate measurements from a high-quality gyroscope. 

The associated error state dynamics arc 

x = Fx + Gu + w (10) 
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Full State Update In an MEKF, the full-state estimate update for the angular rate is made addi- 
tively (u> = Slo ± u>~), and multiplicatively for the attitude quaternion (q + = Sq <g> q~). 4 

Performance Telemetered flight data shows unbiased measurement residuals of roughly 50-arcsec 
(3a) for transverse attitude errors (Af) x y ) and 500-arsec (3a) about the spin-axis. Those numbers 
are consistent with attitude estimation errors of 40 and 1 10-arcsec (3a) respectively (versus sim- 
ulated truth). Similarly, the associated angular rate estimation errors arc 0.01 and 0.03-deg/sec 
(3a) during non-thrusting periods — with errors peaking 3-6 times greater while thrusting. A more 
detailed examination of the MEKF’s performance is planned for a future publication. 

Acceleration Measurement System (AMS) 

The Acceleration Measurement System (AMS), which was manufactured by ZIN Technologies 
of Cleveland, Ohio, provides three-axis acceleration measurements during orbit-adjustments and 
integrates these samples in order to help determine the net velocity change (Av) imparted by these 
maneuvers. The AMS incorporates internally redundant electronics within a single enclosure with 
each side interfacing with its own set of three single-axis orthogonally mounted Honeywell QA3000 
accelerometer sensors. The AMS electronics include high-rate (100 kHz) analog-to-digital signal 
conversion. An embedded processor down-samples and filters the acceleration data. The AMS has 
a dynamic range of greater than ±25,000 pg, a resolution of less than 1 pg, a short-term, la bias 
stability over a twelve hour period of better than 1 pg, and an effective bandwidth of 250 Hz. 

Accelerometer Dynamical Model 

An accelerometer behaves dynamically as if it were a small proof-mass attached to a central-body 
by relatively stiff tri-axial springs. The device then generates a measurement signal proportional to 



a (small) displacement £ relative to its non-accelerating (rest) state. This conceptual model is suffi- 
cient for most applications regardless of the actual device’s internal construction (e.g. cantilevered 
beam, electro-magnetic re-balance loop, etc.). 

Given the accelerometer’s proof-mass position R p with respect to an inertially-fixed origin is 
equivalent to the location of a body-fixed origin’s inertial position (R 0 ) summed with a local vector 
(r ( [) that places the device’s in body-fixed coordinates, plus the proof-mass’ deflection from its 
rest-state. This relationship is expressed as the equation 


R p — R 0 + A b (r d + £) (13) 

r p 

where the symbol A is used generally in this paper to denote an attitude direction-cosine matrix, 
and in this specific instance is the transformation from the body-fixed to inertial frame. Twice 
differentiating Eq. (13) with respect to time produces the kinematic acceleration of the device-mass 
rrip with respect to the inertial frame as 

d 2 

V p = — Rp = V 0 + A (u) x Yp + u) x v.p + 2 • ui x rp + bp) (14) 

(it i^b 

where u> is the angular velocity of the body-fixed frame with respect to the inertial frame, expressed 
in body coordinates. Assuming the tiny proof-mass deflections are much smaller than the vector 
locating the device in the body-frame, then r p ~ r d . If we also neglect structural flexure of the 
base-body, then rigidity implies r d = r d = 0, and the kinematic expression of Eq. (14) reduces 
even further. Substituting for the acceleration at the sensor-head into Newton ’s second law of motion, 
equates the total proof-mass acceleration to the sum of the external forces acting on it as follows 
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where a grav is the acceleration due to a gravitational field. Eq. (16) also moved r p = £ to the 
left-hand side of the equality, and uses Newton ’s third law to replace the proof-mass acceleration 
| with the tri-axial spring-force ( c d is the damping coefficient of the spring, and kd its stiffness). 
If the dynamics-of-interest are far enough below the bandwidth of the sensor (V G , A < ^), the 

accelerometer’s internal response will have damped out (£ — > 0 ), and the device will produce a 
steady-state output 


a d = = a (v o - a grav ) + u x r d + u x u x r d (17) 

Note that the only external force the device’s proof-mass experiences is gravity, because it is a 
field-effect. All other mechanical disturbances to the system must be transmitted through the base- 
body via the tri-axial spring. Prior to launch, resting on a clean-room floor, the accelerometer 
“sees” gravity on its output because the base-structure is at rest (V 0 ~ 0, unless of course there 
is an earthquake or heavy-footed technician nearby). In orbit, the expression for the base-body’s 
acceleration with respect to inertial can be obtained by differentiating the linear momentum P of 
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the base-body with respect to the inertial frame, 


d 2 
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where r c is the position of the base-body’s center-of-mass (CM) expressed in the body-fixed frame, 
m is the spacecraft mass, and 5v rn is a relative exhaust velocity of any expelled mass (e.g. thrusters). 
Even though a rigid-base was assumed in the derivation of the proof-mass motion, choosing not to 
nullify the derivatives of the CM at this point provides a placeholder for multi-body effects such 
as fuel-slosh and appendage motion, so they remain intact for now. Applying Newton ’s second law 
and identifying the mass-expulsion as the body-fixed thruster-force (fj), Eq. (19) now takes the form 

Vo = — y^F b - A (r c + cj x r c + 2cj x r c + u; x cj x r c ) (20) 

= a grav + A ( + — - r c - u x r c - 2 • u> x r c - ^ x a; x r c ) (21) 

i^b \ m m J 

Assuming that the accelerometer proof-mass is much smaller than that of the base -body ( m Y: / <C m), 
the spring reaction-force — £ can be ignored without inducing significant error. Finally, combining 
Eq. (21) with Eq. (17) yields a complete expression for an accelerometer truth-model 

a d = — + d; x (r d - r c ) +lo x uj x (r d - r c ) - (2 • u x r c + r c ) (22) 

m ' v ' 2 ✓ 

r cd multi-body effects 

It its interesting to note the cancellation of the gravitational acceleration terms between the two equa- 
tions. In truth, this only holds if the gravitational field is completely uniform. However, with the 
semi-major axis of the MMS orbit varying from 1.2-12 Earth radii, the so-called gravity-gradient 
effect (defined as the difference in the gravitational field across the spacecraft) is negligible. Specif- 
ically, using those inverse-squared distances with a CM-to-sensor displacement (r cd ) of only one 
meter produces a maximum effect of 0.18 fig at perigee and 2 x HP 4 fig at apogee — sufficiently 
below the threshold of the MMS accelerometer. 

Measurement Equation 

A discrete-sample at time tj. of the acceleration from the n th sensor-head, a n , is modeled as 

Qn(ffc) — (Id - ^ ) • U n a<2 -f- b n Tf n (23) 

where 5 Kb is a scale factor error, u is a unit vector along the true sensing axis of the n th sensor 
head, b n is an intrinsic electro-mechanical bias, and r/ n is sensor noise (possibly non-white). All of 
the introduced quantities (5k, u, b , ?/) arc both time-varying and functions of temperature. A vector 
quantity for the sampled acceleration (a/,) in the body-fixed frame can be reconstructed from the 
n (> 3) measurements using the pseudo-inverse of a user-supplied orthogonality matrix O, which 
ideally cancels the true alignment matrix U = [ui U 2 ... u n ] T and scale-factor errors in a 
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holds true. In the case of MMS, the number of sensor heads is three (n = 3) and — through ex- 
tensive ground calibration and on-orbit temperature control — they are believed to be nearly or- 
thogonal. For the sake of simplicity in the explication that follows, it is assumed the scale factor, 
non-orthogonality, and alignment errors arc essentially zero (although the full effects of the errors 
arc bounded and included in the subsequent Monte Carlo analysis). With the stated assumptions, U 
and O arc both the identity matrix I, and Eq. (24) reduces to 

a k = a d + b + r] (25) 


VELOCITY ESTIMATION 


The quantity of interest from a formation-maintenance perspective is not the acceleration per 
se, but the change in velocity of the spacecraft’s center-of-mass due to thrusting. Analogous to a 
rate-integrating-gyro for attitude dynamics, the AMS’s primary function for the mission is to act 
as an acceleration-integrating-accelerometer during orbital-adjustments. Combining Eq. (25) with 
Eq. (22), transforming to the inertial frame, and integrating over the time-interval (t \ , t 2 ) produces 
the following relationship 



A a k dr 
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(26) 


Evaluating the definite integral, and denoting the integrated thrust as the perturbed velocity of the 
spacecraft’s center-of-mass Sv c (t) the relationship takes the form 




A a k dr = Sv c (t 2 ) - Sv c (ti) + A [u(t 2 )] x r cd (t 2 ) - A [cj(fi)] x r cd (fi) (27) 
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where the attitude matrix subscripts b\ and b 2 refer to the body’s orientation at times t \ and t 2 re- 
spectively. Delaying the expansion of the definite integral limits to compress terms, we arrive at 
an important result. On the left-hand side of the following expression is the true quantity to be 
controlled, Av c , and on the right an expression for it in terms of ideal integrated sensor measure- 
ments and some additional “error” sources that must be removed (estimated/compensated) in order 
to achieve precise maneuvering. 

ft2 r 'l t-2 ft2 ft2 

Av c (ti,f 2 )= / Aa k dr+\ Au x r cd - -4r c > - Ah dr- Arjdr (28) 
v ' J t\ b l ^ b J t, Jt! l *~ b Jt\ l *~ b 
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Methods for estimating the terms of the right-hand side of Eq. (28) will be detailed in the sections 
that follow. 
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Measurement Integral 


The first term on the right-hand side of Eq. (28), the measurement integral, is obtained by sam- 
pling the analog accelerometer (e.g. QA 3000) at a sufficiently high rate to capture all the relevant 
dynamics of the maneuvering spacecraft. In the case of MMS, the AMS electronics samples well 
above 1 kHz, and then applies appropriate decimation, anti-aliasing and noise reduction filters to 
deliver a “clean” 1 kHz digital signal. The discrete AMS output has a bandwidth of 250 Hz, and less 
than ±5 ppm of gain and/or ±1° of phase distortion in the pass-band. Since the MMS ACS oper- 
ates on a 4 Hz control cycle, the measurement integral contains 250 subsamples that are combined 
to generate a single velocity-increment. Outside the AMS, in the ACS flight software, the subin- 
terval sample is transformed into the inertial frame, and summed to produce the total maneuver 
velocity-change estimate for closed-loop control. 

Frame Rotation Compensation Because the summation of the subsamples must occur in a uni- 
form reference frame, each subsample is transformed forward in time to what is assumed will be 
the frame of the final (k = 250) sample. This transformation is based upon the current estimate 
of the spacecraft’s expected attitude motion. The propagation assumes the base-body’s angular rate 
estimate u>(/i) is approximately constant over the full quarter-second sample interval. In the par- 
lance of strap-down inertial navigation systems, this type of compensation to the velocity update is 
known as a sculling correction. 6 ’ 7 (For reference, the dual term applied to gyro measurements for 
attitude updates is referred to as a coning correction.) 


Recalling the measurement integral term of Eq. (28), its discrete approximation is 

t2 250 

/ A (T)a k dT& 4 V 4 (u>,k) - a* • (t k - 

J £ i<—b 250 z J 250 <— k 


(29) 


k = 1 


where the attitude matrix update is obtained from an Euler axis/angle parametrization 5 


A (e, <hfc) = I — sin 4>fce x + (1 — cos <I>fc) e x e x (30) 

250<-fc 


where the Euler axis e = u>/ ||a> || is a unit vector in the direction of the angular rate u expressed 
in body-fixed coordinates, and the Euler angle $ is a function of the time-difference between the 
current (k th ) and final (250 ,/l ) sample, 


$ k = H&ll • A t k 




250 -k 
1000 


(31) 


Recall that u> is the axis of rotation and assumed constant over the sub-interval. It is therefore the 
same in all (instantaneous) subsample coordinate frames. 


Centripetal Compensation 


An estimate for the kinematic effect that results from having the accelerometer mounted with an 
offset from the structure’s spin-center (i.e. the center-of-mass) is obtained by replacing the values 
in the dynamical terms with current-best-estimates. In the case of the attitude transforms A and the 
angular rate u>, these estimates come directly from the on-board attitude estimate whose accuracy is 
governed by the star sensor measurements and MEKF solution. The only other estimate needed is 
the offset of the accelerometer sensor heads from the system's center-of-mass (r 
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(32) 



It is worth emphasizing the efficiency of this result. Since the centripetal correction is applied to the 
integrated acceleration — and is itself a velocity-correction — it is applicable over any timespan with 
just this simple two-term computation. Even if the time-interval (ti , ^ 2 ) spans multiple revolutions 
of the base-body, this correction is exact to within the accuracy of our estimates. 

Center-of-Mass Knowledge Knowledge of the location of the center-of-mass of a 939-kg ob- 
servatory with eight appendages and slightly over 41 1-kg of hydrazine fuel proved to be one of the 
most challenging aspects of the mission. The risk to maneuvering accuracy was mitigated by a com- 
bination of techniques. First, the mechanical book-keeping on the structure was accomplished to an 
accuracy impractical prior to the extensive use of solid CAD models of an assembly. 8 Secondly, the 
calculated mass properties for the fully-stowed (dry) spacecraft were verified prior to launch using a 
spin-balance facility Miller table. Additionally, a number of propellant loading/unloading tests and 
Monte Carlo simulations were performed to bound the final range of “shapes” the fuel — constrained 
by an elastomer-diaphragm — might take in the four cylindrical tanks prior to launch. Lastly, a se- 
ries of on-orbit calibration maneuvers were performed while capturing the 1 kHz accelerometer 
data. The telemetered data was processed on the ground using an augmented MEKF for parameter 
estimation/system-identification. Ultimately, it is believed the CM statistical knowledge-accuracy 
achieved was better than 1 -millimeter (la) for the observatories prior to Axial Double Probe and 
Spin-plane Double Probe (SDP) deployment. It can be seen by inspection of Eq. (22) that if ele- 
ments of the angular rate are constant, a fixed error in center-of-mass knowledge Sr c would yield an 
error u x u x r d that is indistinguishable from a sensor bias. A future publication is being prepared 
that will explore this and the larger topic of MMS system identification in greater depth. 

Sensor Head Offsets A subtlety exists with centripetal compensation that is due to the multiple 
sensor-heads (n > 3) that are required for a full vector solution. It is physically impossible for the 
uni-axial heads to be exactly co-located, so the centripetal compensation must account for each of 
the small head-offsets. MMS was fortunate that this correction could be applied after integration — to 
the accumulated velocity-output. In general, when the sensing axes arc not all mutually orthogo- 
nal, the head-offset correction would need to be applied to the acceleration measurements directly, 
before they arc collapsed in to a single vector-quantity via Eq. (24). 


Multi-body Effects 

Multi-body effects — more precisely the internal flexure of the deployable appendages and flow 
of fuel — defeats most attempts at simple analytic approximation in the sense of the classical rotat- 
ing rigid-body equations. In the course of analyzing MMS, a number of linearized motion models 
were used to gain insight regarding the system modes, 9 but ultimately it was the non-linear time- 
domain simulations of the spacecraft that were used to predict the multi-body response of the sys- 
tem. Specifically, the boom- appendage and fuel-slosh dynamics were modeled using multiple rigid 
bodies connected by rotational joints. 10 

However, regardless of the complexity of the flexible -body dynamic motion, it is true that its 
effect on the spacecraft is non-propulsive. As a result, the expected value of the center-of-mass 
motion term in Eq. (28) will — over a sufficiently long interval — integrate to zero. With regards to 
acceleration feedback correction/compensation, CM-motion during a maneuver was ignored in the 
MMS estimator-design with the caveat that the nutation and wire -boom deflections were minimized 
during the maneuver using interleaved angular momentum control-pulses for which the moniker 
DVDH was coined. 
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One area in which the boom flexibility appears directly is the MEKF. It was shown analytically 
and via simulation that the state-propagation is improved if an “effective inertia” is used for the 
filter’s rigid-body approximation of the rotational dynamics. Following the method of Zimbelman 
and Walker, 1 1 the effective inertia about the spin-axis I e ff used in the MEKF is defined as 

h ft = /base + ( /sdp (33) 

v r sdp + LsdpJ 

where r sc | p is the radial distance from the spin-axis to the SDP pivot-point, L sdp is the distance along 
the wire -boom to the SDP’s center-of-mass, and I sc j p is the SDP’s second mass-moment of inertia 
about its own center-of-mass as if it were a rigid rod. A similar expression also is used for the 
effective inertia’s about the spacecraft’s transverse axes. 

Bias Estimation 

A dominant factor in the overall maneuver execution error of the system is the integral of the 
intrinsic electro-mechanical biases. The MMS sensor has the potential for a non-stationary bias on 
the order of ±2000 pg, with a one-year composite bias repeatability of 105 fig . 1 2 However, over 
a much shorter period of twelve hours, the bias stability (in a temperature controlled environment) 
is less than ±1 pg. As a result, over the duration of a maneuver (< 90 minutes for MMS), the 
bias is nearly constant and is compensated internally in the AMS by subtracting a per-maneuver 
bias estimate b from the acceleration measurement subsamples before they are added to the running 
velocity total. The bias-residual is defined as the error in the bias estimate (5b = b — b). The 
estimate b is itself obtained from an independent filter internal to the AMS. Just like an attitude 
Kalman filter uses the inertial reference of star sensor measurements to estimate the biases in gyros, 
the inertial reference of a non-thrusting period of spacecraft orbit (Av c = 0 ) is sufficient to deter- 
mine the intrinsic bias of the accelerometers. In the case of MMS, a 6 ,/ ' -order, cascaded (bi-quad) 
low-pass filter was used with a frequency break-point set at 0.004 Hz (4.2 minutes). Due to its long 
time-constant, the bias-residual filter requires 30-minutes of quiescent data prior to each maneuver 
in order to settle to within ±1 pg of its steady-state solution. In order to speed-up convergence and 
also compensate for any base-body nutation, the expected value (or mean) of the gyro-dynamic bias 
is subtracted from the bias-residual filter’s input. The mean gyro-dynamic bias is determined by 
taking the expected value of Eq. (25) and expanding for a,; using Eq. (22). 

E [a fc ] = E [<± x ] r cd + E [u> x u> x ] r cd - E [2 • u> x r c + r c ] + E [b] + E fa] (34) 

Recognizing that the expected value of the zero-mean noise and multi-body (internal flexure) terms 
are both necessarily zero, E [b] is the bias estimate b being sought, and add to that E [a/.] is ap- 
proximately the low-pass bias-residual filter’s output, to obtain the following new expression 

b = E [afc] - A ( E[u £] A r cd + E [c^>„ x ] A r cd ) (35) 

b<—p y p^b p^b J 

bias filter output 

where u> p is the angular rate with respect to inertial of a body-fixed frame p aligned with the space- 
craft’s principal axes. The principal-axis frame was introduced here because, for an axisymmetric 
major axis spinner, there is a well known analytic solution for its motion 13 

u\(t) = ttq(0) cos Xt — 0 ^ 2 (0) sin At (36) 

L 02 (t) = wi(0) sin Xt + u) 2 (0) cos Xt (37) 

W 3 (*) = w 3 (0) (38) 
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where A is the nutation frequency derived from the principal moments of inertia (/■; > I 2 > I\) and 
(constant) major-axis spin-rate as follows 


A = (39) 

h 

where the transverse moment of inertia I t = I\ = I 2 in an axisymmetric body (or can be approxi- 
mated using the geometric mean 7/ ~ y 7 1\ • I 2 in the nearly axisymmetric case). From the coupled 
motion equations, it can be determined by inspection that E \cb p ] = 0 and 


E [^p ] A r cd 

p^b 


-w 2 s ■ ri + E [wiwa] r 2 - E [uj 2 ] n 
-w 2 ■ r 2 - E [ojf\ r 2 - E [u > iuj 2 ] r\ 
~ ( E [A] + E [^i] ) r 3 


w 2 ■ r\ 


r 19 1 

-2 u t - ri 

w 2 ■ r 2 

+ 

1 2 

~2 U t ~ r 2 

0 


1 

1 

£ 

^ to 

-i 

CO 

1 


spin bias nutation bias 


(40) 


(41) 


where uj s = oj:> ( 0 ) and uy = \J wf + w \ arc the magnitude of the spin and transverse rates in the 
principal axis frame respectively. The first term in Eq. (41) for MMS has a magnitude of 3950 
fig, however the nutation term is small — on the order of a few micro-g — and to the best of our 
knowledge, had not been specifically identified in the published inertial-navigation literature. 

As an additional consideration, signal-aliasing of the wide -band noise must be carefully managed 
to prevent corruption of the final bias estimate. For this reason, the bias residual filter was located 
inside the AMS, leveraging its ability to process the high-frequency data directly. 

Principal Axis Alignment Since the average instantaneous spin of the spacecraft is about the 
major principal axis, the orientation of the principal-axes frame with respect to the body-fixed axes 
is estimated prior to each maneuver using a low-pass filter. This estimator resides in the ACS flight- 
software and operates on the body’s rate-estimate Cb produced by the attitude determination MEKF 
at 4 Hz. Recalling Poinsot’s construction 5 - 13 for an axisymmetric body, the body-cone angle (3 can 
be determined from the components of the filtered angular rate uij as 


(3 = tan 






(42) 


This angle is the second rotation in a classic 3-1-3 ((/> — (3 — ip) Euler angle sequence that transforms 
from the body to the principal axes frame. The first rotation in the sequence, ((>, is 

<t> = tan -1 ^ (43) 

\Uf*J 2 

and the final clocking angle w of the minor and intermediate axes about the major axis may be 
chosen arbitrarily, since it does not affect the bias estimation calculation in Eq. (41) (nor do the 
exact magnitudes of the principal moments). The final set of rotations can then be used to construct 
the direction cosine matrix A or equivalent quaternion. 

p<—b 
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Effects of Bias Integration in a Rotating Frame An important consideration for maneuver knowl- 
edge error is what happens on a spinning platform to a bias that is constant in the body-frame after 
it is transformed to the inertial-frame and integrated. For a pure z-axis rotation, we have 


nt+At 

v = / A (r) ■ b • dr 

Jt 


rt+At 


COS UJ Z T 
— sin Ld z T 
0 


sin oj z t 
cos U 1 Z T 

0 



sin ui z At + ^ (1 — cosu^A t) 
(cos uj z At — 1) — ^ sin u> z A t 
6, At 


O' 


t>X 

0 


by 

1 


b z _ 


dr 


(44) 


(45) 


This simple result shows that the effect of integrating bias residuals along the two transverse axes 
is both cyclical and bounded. Flowever, the axial Av knowledge error due to integrated bias es- 
timation error grows linearly with time. In contrast, the RMS of noise-driven velocity random 
walk (discussed next) only increases as the \fK t. For all practical MMS maneuver durations, the 
integrated axial bias residual is the dominant source of error. This qualitative difference due to 
direction (axial/radial) is revealed in Figure 3 during a portion of a high-fidelity simulation of the 
fully-deployed spacecraft executing a 5 m/sec maneuver. 



Figure 3. Simulation of Axial and Radial Knowledge Errors 


In order for MMS to meet its mission requirements, it was expected that the composite bias- 
I'csidual/CM -knowledge-bias error would be < 1 /ifj (la) prior to a maneuver. 


Velocity Random Walk 


The final right-hand term in Eq. (28) is the result of conceptually integrating the wide -band sta- 
tionary noise of the accelerometer. This error source is assumed to be driven by an unspecified 
electrical-thermal source. An estimate of its effect on velocity can be determined by taking it ex- 
pected value 


E[ Avf = E 


rt2 


A ■ T]dT 



E[rj\ ■ dr = 0 


(46) 
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where the body to inertial transformation A may be dropped because the noise is temporally uncor- 
related. While the expected value of the noise is zero-mean — and will not bias the estimate of the 
velocity — the variance (mean-squared value) of the integrated noise is 

rt2 -1 roc 

E [ Av l} ~ / E [v 2 ] dr = — S(juj) du> ■ dr (47) 

Jtl Jt 1 Z7T 7-00 

where S(ju ) is the power spectral density (PSD) of the noise, and a; is a frequency in rad/s — not 
to be confused with the base-body angular rate u. If the input is assumed to be band-limited white 
noise, then the PSD is a constant A/2 within the pass-band (/ hz )), and zero otherwise. The velocity 
variance becomes 

/•At /-At i /•27r/ hz a rAt 

E[ Av%] = / E [r? 2 ] dr = — — du; ■ dr = A f hz dr = A / hz At (48) 

Jo Jo 77T 7-27r/ hz 7 Jo 


The QA3000 specification asserts the root-mean-squared noise (a rms = yj E\ Av^]) should be less 
than 8 //// rms in the frequency range 0-10 Hz, and 80 /t/y rms from 10-500 Hz. This implies that if 
a (single-sided) periodogram of acceleration-samples over a large time-interval was computed, its 
expected magnitude {A) over a given frequency range (e.g. A/h z = 10 Hz) should be less than 


4lo-10hz < 


(Orms) 2 

A/hz 


sW =64 /V 

10 Hz ' Hz 


615.9 


Hz 


(49) 


Figure 4 shows this is indeed the case for a sample taken from an non-maneuvering MMS-1. With 



frequency [Hz] 


Figure 4. PSD of MMS-1 AMS 1 kHz Acceleration Data 

this base-band noise input driving an integrator-plant, the velocity random walk is then 

(Ar^) rms < yj Ao-tohz • Af = (24.8 ^ • VHz) • VAt (50) 

It can be shown that the noise-response is dominated by the content in the base-band (0-10 Hz) due 
to low-pass nature of the integrator (-20 dB/decade). This mitigation strategy for this error-source 
is to limit size of the velocity-accumulation window tightly on the thrusting portion of a maneuver. 


Other Error Sources 

Non-Linearity At operational temperatures, the non-linearity error (of the scale-factor) over the 
QA3000's dynamic range was calculated to be less than 150 ppm. 
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Thermal Effects. Nearly all elements of the AMS exhibit some degree of temperature depen- 
dence that could affect performance. For MMS, these error were minimized by actively controlling 
the thermal environment — limiting variations during maneuvers to with ± 1-degree Celsius. Addi- 
tionally, a number of temperature corrections arc applied internally to the acceleration samples by 
the AMS based its own internal thermistors with 0.1° resolution. 

Scale Factor Through careful component selection, and by performing extensive ground cali- 
bration of the AMS in an isolated environment (NASA Glenn Research Center, Plum Brook Facil- 
ity), ZIN was able to limit expected scaling factor errors to less than 872 ppm. 

Non-Orthogonality / Sensor Alignment Alignment of the individual sensor-head axes was con- 
trolled by machine tolerances, thermal compensation during flight, and verified by ground metrol- 
ogy within an accuracy of 1 15-arcsec per axis. This error — contributing less than 0.2-ppm of scale- 
factor error and mapping directly to directional error (115\/3 arcsec = 0.055°) — does not dominate. 

Quantization The resolution of the sensor is better than 0.8 fig. As a result, quantization through 
discretization is not considered a significant source of error. 

Frequency Response / Structural Dynamics The AMS is mounted on a shock/vibration isolation 
system designed by MOOG/CSA Engineering. Nevertheless, the structural response seen at the 
accelerometers due to thruster firings is very pronounced — at times, greater than 75% of saturation. 
The MMS structure exhibits a number of modes in the 30-300 Hz range that pollute the AMS 
acceleration output. Fortunately, due the low-pass characteristics of integration, the ripple-effect 
on maneuver velocity knowledge is below 1 mm/sec (see Figure 5), and the small error can be 
absorbed in the performance budget. Furthermore, the structural-response damps quickly — only 
three-seconds of padding is required at the end of each maneuver to account for ring-down. 



783.25 783.5 783.75 784 784.25 784.5 784.75 785 783.5 784 784.5 785 

Time [sec] Time [sec] 

Figure 5. Structural Modes Effect on Accumulated Velocity (telemetry) 

Summary of Velocity Estimation Method 

Having first briefly explained how the spacecraft constructs its attitude and rate knowledge from 
star sensor measurements, details were then provided on an acceleration measurement model. It was 
shown that directly integrating the acceleration-samples would not produce an accurate estimate of 
the spacecraft’s center-of-mass velocity change, so the necessary corrections were derived in detail. 
All of the compensations explained may be categorized as follows: 

• bias estimation (pre-maneuver) — The ACS and AMS use the angular rate with long duration 
low-pass filters to arrive at an accurate intrinsic sensor-head bias estimate. 

• bias and sculling correction — The AMS uses the bias and angular rate estimates to correct 
and sum the (1 kHz) acceleration-samples in order to produce a (4 Hz) incremental velocity. 
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• centripetal correction and inertial summation-The angular rate is used to apply centripetal 
correction to the (4 Hz) incremental velocity output, and the total inertial (ECI J2000) velocity- 
change imparted to the spacecraft is accumulated. 

ATTITUDE, ANGULAR MOMENTUM, AND TRANSLATIONAL CONTROL 

Actuators (thrusters) 

Both attitude and orbital control of the observatories arc accomplished using twelve hydrazine 
mono-propellant thrusters — four AMPAC 1-lbf (4.4 N) thrusters arc directed axially (±z) and eight 
Aerojet 4-lbf (17.8 N) thrusters radially (±y). The minimum impulse bit for the Aerojet designed 
thrusters ranges throughout the mission from 0.13-0.26 N-m-sec, which corresponds to a 20 mil- 
lisecond pulse. In addition to a continuous firing mode, the engine-value-driver is capable of mod- 
ulating the duty-cycle of the output on 8 Hz boundaries. This pulse-width-modulation allows for 
balancing of torques during translational maneuvers that minimizes attitude perturbations. 

Velocity Control Vector Definitions 

With the definition for the truth-states of Eq. (28) in-hand, a number of other important control- 
vector quantities may also be derived. First, the imparted velocity estimate Av c is defined as the 
expected value of the true velocity-change, 

Av c (fi,£ 2 ) = E[Av c (t 1 ,t 2 )} (51) 

The maneuver knowledge-error (<5v) is the difference between the estimate and truth 

8v(h, t 2 ) = Av c (fi, t 2 ) - Av c (fi, t 2 ) (52) 

The true performance metric for the system is the control-error 8v, which is the difference between 
the commanded velocity-change (target) and the truth 

8v(ti, t 2 ) = Av c (ti,i 2 ) - Av tgt (ti, t 2 ) (53) 

Finally, the definition of the servo-error (5v servo ) is the difference between the target and estimate, 

<5v serV o = Av tgt (/i,f 2 ) - Av c (ti,i 2 ) (54) 

and is the quantity that the controller actively regulates. 

Velocity Control 

The AV controller follows the design principals of a classic tracker — with both a time-varying 
target, and velocity-estimate feedback. A predetermined velocity-change profile is uploaded to the 
spacecraft prior to each maneuver. The trajectory is in the form of a piece-wise linear look-up 
table — AV in ECI J2000 versus spacecraft time. There arc very few constraints on the nature of the 
profile. For example, it may be non-mono tonic ally increasing in magnitude, and/or it may possess 
arbitrary spacial curvature. However, it must be constructed with a care towards the saturation limits 
of the thrusters or the system will perpetually lag behind the target — potentially failing to regulate 
the servo-error sufficiently within a given maneuver window. The controller is entirely responsible 
for the mapping of velocity-commands — specified in the inertial frame — to the appropriate group- 
ing of thrusters (i.e. “bank”) capable of delivering an ideally torque-less translational force to the 
spacecraft. The projection of the servo-error into cylindrical coordinates determines if an axial or 
radial bank is a proper candidate for actuation. 
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Axial Thrusting On a spinning spacecraft, axially-aligned thrust is conceptually simpler. Since 
the spacecraft’s center-of-mass radial offset from the spin-axis is nominally zero, no torque-balancing 
via pulse-modulation is attempted. However, since the axial thrust is always well-aligned with the 
axial projection of servo-error, in order to reduce the amount of valve-cycling and chattering be- 
tween the top and bottom (±z) banks, the controller is designed to utilize a dead-hand, and limits 
commanding to roughly once every quarter spin period. Since translational and rotational control is 
interleaved (and translational-control is given priority), this “packeting” of axial-actuation also of- 
fers a potential benefit by freeing those thrusters for full momentum/attitude-control. The qualifier 
“full” was applied to the previous statement because axial off-pulsing of a single thruster in a pair 
is also sometimes necessary — especially when a AV profile contains a dominant axial-component 
and/or the (diaphragm-constrained) fuel-mass is unbalanced about the spin-axis. 

Radial Thrusting The MMS radial thrusters 
are arranged in two banks of four, each di- 
rected along the positive and negative y-axis 
respectively. Pulse-Width-Modulation (on 8 
Hz boundaries) of the upper or lower pairs in 
each quad facilitates torque-balancing. This 
arrangement offers two opportunities per spin 
for thrust along a given inertially-fixed direc- 
tion (roughly every 10-seconds). In order to 
avoid inducing directional error, the radial ac- 
tuation is centered around the servo-error’s pro- 
jection in the spin-plane. Because the controller 
tracks a time-varying velocity-profile, the burn- 
centering problem is non-trivial. The precise 
timing affects the burn-vector and vice versa — 
therefore an iteration solver is employed to re- 
solve target-centering down to a millisecond resolution. Explicitly, for a given servo-error, the 
current rate estimate can be used to predict the time at which the thrusters will align with the center 
of the required burn-arc. However, at the predicted center-time, the target will have most likely 
shifted — resulting in a slightly different servo-error projection, ad infinitum. Fortunately, conver- 
gence is typically achieved within a few iterations. Also, to prevent the controller from constantly 
lagging it’s target, a quarter spin-period “look-ahead” strategy is employed. As shown in Figure 
6 this 5 -seconds of time-bias keeps the servo-error straddling the velocity -profile, and speeds up 
entry into a fine-tuning ’’trim” phase at the end of each maneuver. Similar the axial control strategy, 
dead-hands arc employed to reduce the wasting of fuel due to chatter. 

Momentum Control The spin-axis inertial pointing direction is important to the MMS science 
instrument suite. For this reason, each observatory performs a momentum-adjustment (AH) imme- 
diately following every AV. While the specifics of the AH control-logic are not examined here (see 
the reference 14 ), it is also important to hold pointing, minimize nutation, and maintain spin during a 
AV-maneuver — not only for the instruments, but also for accurate velocity-estimation. Momentum 
control during a AV maneuver (DVDH) employs the same Fyapunov-based controller as the stand- 
alone AH mode, with two main differences — the aforementioned off-pulsing of axial-thrusters, and 
the frequency of the control (e.g. while DVDH is only allowed to actuate up to once every 3-seconds 
in order to avoid structural resonance with the two ADP booms, AH-mode is typically slower). 



Figure 6. Velocity Controller Telemetry (FI-116) 


16 


SYSTEM ROBUSTNESS 


Two major concerns of control system design are robustness in stability and performance in the 
presence of plant uncertainty. For MMS, exhaustive Monte Carlo simulations were used to address 
these twin concerns. Following the statistical methodology of Flanson and Beard, 15 99% confidence 
criteria (1% consumer risk) was selected that requires zero failures over a sample-size of 3410 runs. 

This criteria was tested repeatedly using GSFC’s Freespace Simulation Environment 16 to sta- 
tistically vary over 250 model parameters — resulting in hundreds of thousands of time-domain 
simulations of maneuvers at full model-fidelity. Figure 7 depicts one example of these results, 
and is annotated with the performance criteria used to determine execution error acceptance for a 
formation-maintenance class of closed-loop maneuver (0-10 m/sec). The scatter plots of maneu- 
ver magnitude and directional error from 3500 runs shows ample margin versus the 3cr closed-loop 
(AMS) requirement (in cyan). 
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Figure 7. Monte Carlo Results for Formation Maintenance Maneuvers 
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MMS MISSION PERFORMANCE 


Open-loop Maneuver 

During early observatory commissioning, a series of perigee -raising maneuvers were performed 
to lift the orbits above the launch vehicle insertion-altitude and reduce atmospheric drag. Because 
these maneuvers were relatively large, and occurred before the AMS performance had been fully 
verified, they were executed open-loop (i.e. without accelerometer feedback). Although it was 
not in the control-loop, the AMS was still available to monitor the performance of the open-loop 
controller. The data shows that, despite a qualification test-program that attempted to characterize 
thruster-performance prior to launch, the propellant-system still exceeded its expected-value by 
approximately 8% for the radials (nominally 4-lbf). This only underlies how difficult it is to predict 
total system behavior during mission development, and emphasizes the need for closed-loop control 
when a high degree of accuracy is required. One example of a difficult problem to characterize is 
performance during warm-up of the hydrazine catalyst-beds. Flight and ground-test data showed 
that, during the first 20-seconds after a cold-start, thrust could start as much as 15% lower than the 
final steady-state output. The plot of flight-telemetry shown in Figure 8 emphasizes the divergence 
between the open and closed-loop Av c estimates for an MMS-1 perigee -raise maneuver. 
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Figure 8. MMS-1 Open-loop Maneuver Performance Comparison to Closed-loop Estimate 
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Closed-loop Maneuvers 

The best avenue for determining the true closed-loop performance is through post-maneuver orbit 
determination. The GSFC developed NAVIGATOR 17 Global Positioning System (GPS) and GPS- 
Enhanced Onboard Navigation System (GEONS) units arc collectively providing autonomous orbit 
determination with a record-breaking 5 meters (3er) semi-major axis (SMA) accuracy. 18 After the 
first few closed-loop formation initialization maneuvers for MMS commissioning, the following 
performance has been verified 


Maneuver 

(DOY) 

Obs 

ID 

Final Target 
Magnitude 

GEONS Solution 
Semi-major Axis 

Final 

Servo-Error 

AMS Bias Estimate 

Gg) 

mm/s 

A-error 

mm/s 

% target 

X 

Y 

Z 


1 

118.6 

-1.14% 

1.5 

1 .25% 

114.7 

78.9 

49.6 

GS-095 

2 

18.3 

-0.57% 

1.0 

5.73% 

94.3 

93.9 

47.3 

(163) 

3 

46.9 

-0.73% 

1.1 

2.27% 

75.2 

92.5 

140.1 


4 

77.0 

0.55% 

1.1 

1 .44% 

108.3 

96.1 

125.1 


1 

0 

— 

— 

— 

115.3 

77.4 

49.7 

FI-116 

2 

4077.5 

-0.79% 

1.0 

0.03% 

95.0 

94.0 

47.5 

(184) 

3 

9175.6 

-0.26% 

0.2 

0.00% 

76.9 

94.3 

140.9 


4 

4452.1 

-0.26% 

1.2 

0.03% 

107.2 

93.9 

125.4 

FI-118 

1 

2 

0 

3511.6 

-0.61% 

0.8 

0.02% 

93.7 

94.0 

47.6 

(190) 

3 

4149.7 

-0.18% 

1.3 

0.03% 

76.9 

94.7 

140.8 


4 

6068.7 

-0.27% 

1.3 

0.02% 

106.9 

95.5 

125.3 


While the error in targeted semi-major axis does not equate with controller-error (<5v), SMA is 
directly linked to the orbital speed (e.g. the vis-viva equation). These early results indicate that the 
fleet is maneuvering successfully within the 3cr total mission requirements that arc an amalgamation 
of ACS execution and flight-dynamics planning errors. These results arc particularly encouraging 
considering the sizes of the GS maneuvers, since high percent-accuracy is more challenging for 
small maneuvers. The sample mean [i of these maneuvers is -0.068%, with a sample standard 
deviation a of 0.637%. Based on this relatively small sample size of 10 cases, we can state with a 
90% confidence that the true standard deviation lies in the range of 0.465-1.048%, with the caveat 
that the distribution of maneuver errors is assumed Gaussian. Ultimately, the criteria for a successful 
MMS maneuvering-system design will be the length of time for which a high-quality formation can 
be preserved — which will be measured in the frequency of corrective maneuvers. This will be 
especially evident when the formation moves from its current average separation of approximately 
160 km down to as little as 4 km. 
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